cd "${mystart}/Simulated data\Calibration/Benchmark Data/"

********************************************************************************
******* Compare H, R, du, dp, s to the data *****
********************************************************************************	

* Figure 4
* R profiles
qui{
* looking at R (not by education)
clear
use "Data files/Benchmark.dta"
replace R_recorded=R_recorded-1
replace R=R-1
collapse (mean)   R R_recorded , by(age) fast
		replace R_recorded=R_recorded*100
		replace R=R*100

merge age using "${mystart}\MEPS\MEPS Data Files\Intermediate Data Files\R_ages_ALL.dta"
sort age
replace R3=R3*100
replace R5=R5*100

drop if age>=65
#delimit ;
line R age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) 
|| line  R_recorded age , lwidth(medthick) lcolor(cranberry)   lpattern(dash)
|| line  R5 age , lwidth(medthick) lcolor(edkblue)  lpattern(dash)
  ||, 
  ylabel(0.0(10)50)
  xlabel(25(10)65)
  ytitle("Percent")
xtitle("Age")
   title("High R", color(black))
  legend(label(1 "R Actual, Model") label(2 "R Recorded, Model") label(3 "R Recorded, Data")  )
  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr
graph export "${out_figures}\R_Profile_ALL.eps", replace
}

* Figure 5
* H transitions by education and age, Good to Good, some college
qui{
clear
use "Data files/Benchmark.dta"
sort ID age
gen H_trans=1 if (H==3 ) & H[_n+1]==3  & ID==ID[_n+1]   // stay in good  
replace H_trans=0 if H==3 & (H[_n+1]==2 | H[_n+1]==1)   & ID==ID[_n+1]  // transition to worse health


keep if age>=25 & age<75
gen d_shock=0 if du_record==0 & dp_record==0 & age<65
replace d_shock=1 if (du_record==1 | dp_record==1) & age<65
replace d_shock=0 if du==0 & dp==0 & age>64
replace d_shock=1 if (du==1 | dp==1) & age>64


collapse (mean) GG_model=H_trans, by(education  age_group d_shock) 
rename age_group age
rename education educ_group
sort educ_group age d_shock


merge educ_group age d_shock using "${mystart}\MEPS\MEPS Data Files\Intermediate Data Files\Health Trans Moments data shocks.dta"

drop _merge
reshape wide GG_model GG, i(age d_shock) j(educ_group)

sort age


#delimit ;
line GG_model2  age if d_shock==0, yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) 
|| line GG_model2  age if d_shock==1, yaxis(1) xaxis(1) lwidth(medthick) lcolor(navy) 
|| line GG2  age if d_shock==0, lwidth(medthick) lcolor(cranberry)  lwidth(medthick) lpattern(dash)
|| line GG2 age if d_shock==1, yaxis(1) xaxis(1) lwidth(medthick) lcolor(navy) lpattern(dash)
  ||, 
  ylabel(0.0(.2)1)
  xlabel(25(5)74)
  ytitle("Fraction")
xtitle("Age")
   title("H Transitions from Good to Good, Some College", color(black))
  legend(label(1 "No Recorded d{superscript:p} or d{superscript:u}, Model") label(2 "Recorded d{superscript:p} or d{superscript:u}, Model")   label(3 "No Recorded d{superscript:p} or d{superscript:u}, Data") label(4 "Recorded d{superscript:p} or d{superscript:u}, Data") )
  legend(col(2) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr	
graph export "${out_figures}\Health_Trans_Profile_GG_shock2.eps", replace


}

* Figure 6
*h transitions by education and age
qui{
clear
use "Data files/Benchmark.dta"
sort ID age
gen H_trans=0 if (H==3 ) & H[_n+1]==3  & ID==ID[_n+1]   // stay in good  
replace H_trans=1 if H==3 & (H[_n+1]==2 | H[_n+1]==1)   & ID==ID[_n+1]  // transition to worse health

// from bad to good
gen H_trans1=0 if  (H==2 | H==1) &  (H[_n+1]==2 | H[_n+1]==1)  & ID==ID[_n+1]
replace H_trans1=1 if (H==2 | H==1) &  H[_n+1]==3   & ID==ID[_n+1]

// from poor to better
gen H_trans2=0 if  ( H==1) &  ( H[_n+1]==1)  & ID==ID[_n+1]
replace H_trans2=1 if ( H==1) &  (H[_n+1]==2 |  H[_n+1]==3)  & ID==ID[_n+1]


// from poor to poor
gen H_trans3=0 if  ( H==1) &  ( H[_n+1]!=1)  & ID==ID[_n+1]
replace H_trans3=1 if ( H==1) &  (H[_n+1]==1)  & ID==ID[_n+1]

// from fair to fair
gen H_trans4=0 if  ( H==2) &  ( H[_n+1]!=2)  & ID==ID[_n+1]
replace H_trans4=1 if ( H==2) &  (H[_n+1]==2)  & ID==ID[_n+1]

// from poor to poor
gen H_trans5=0 if  ( H==3) &  ( H[_n+1]!=3)  & ID==ID[_n+1]
replace H_trans5=1 if ( H==3) &  (H[_n+1]==3)  & ID==ID[_n+1]


// from fair/good  to fair/good
gen H_trans6=0 if  ( H==3 | H==2) &  ( H[_n+1]==1)  & ID==ID[_n+1]
replace H_trans6=1 if  ( H==3 | H==2) &  ( H[_n+1]!=1)  & ID==ID[_n+1]


keep if age>=25 & age<75

//drop age_group
//rename age age_group

collapse (mean) GB_model=H_trans BG_model=H_trans1  PG_model=H_trans2 PP_model=H_trans3  FF_model=H_trans4  GG_model=H_trans5 FGP_model=H_trans6, by(education  age_group) 
reshape wide GB_model BG_model PG_model PP_model FF_model GG_model FGP_model, i(age_group) j(education)
rename age_group age
sort age


merge age using "${mystart}\MEPS\MEPS Data Files\Intermediate Data Files\Health Trans Moments data.dta"

sort age
#delimit ;
line GB_model1 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) 
|| line GB_model2 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(navy) 
|| line GB_model3 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(green) 
|| line GB1 age , lwidth(medthick) lcolor(cranberry)  lwidth(medthick) lpattern(dash)
|| line GB2  age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(navy) lpattern(dash)
|| line GB3  age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(green) lpattern(dash)
  ||, 
  ylabel(0.0(.05).4)
  xlabel(25(5)74)
  ytitle("Fraction")
xtitle("Age")
   title("H Transitions from Good to Fair/Poor", color(black))
  legend(label(1 "HS or Less, Model") label(2 "Some College, Model") label(3 "College, Model")  label(4 "HS or Less, Data") label(5 "Some College, Data") label(6 "College, Data"))
  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr	
graph export "${out_figures}\Health_Trans_Profile_GB.eps", replace

sort age
#delimit ;
line BG_model1 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) 
|| line BG_model2 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(navy) 
|| line BG_model3 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(green) 
|| line BG1 age , lwidth(medthick) lcolor(cranberry)  lwidth(medthick) lpattern(dash)
|| line BG2  age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(navy) lpattern(dash)
|| line BG3  age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(green) lpattern(dash)
  ||, 
  ylabel(0.0(.1).7)
  xlabel(25(5)74)
  ytitle("Fraction")
xtitle("Age")
   title("H Transitions from Fair/Poor to Good", color(black))
  legend(label(1 "HS or Less, Model") label(2 "Some College, Model") label(3 "College, Model")  label(4 "HS or Less, Data") label(5 "Some College, Data") label(6 "College, Data"))
  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr	
graph export "${out_figures}\Health_Trans_Profile_BG.eps", replace


sort age
#delimit ;
line PG_model1 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) 
|| line PG_model2 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(navy) 
|| line PG_model3 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(green) 
|| line PG1 age , lwidth(medthick) lcolor(cranberry)  lwidth(medthick) lpattern(dash)
|| line PG2  age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(navy) lpattern(dash)
|| line PG3  age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(green) lpattern(dash)
  ||, 
  ylabel(0.0(.1)1)
  xlabel(25(5)74)
  ytitle("Fraction")
xtitle("Age")
   title("H Transitions from Poor to Fair/Good", color(black))
  legend(label(1 "HS or Less, Model") label(2 "Some College, Model") label(3 "College, Model")  label(4 "HS or Less, Data") label(5 "Some College, Data") label(6 "College, Data"))
  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr	
graph export "${out_figures}\Health_Trans_Profile_PG.eps", replace

gen FGPP_model1 = 1- FGP_model1
gen FGPP_model2 = 1- FGP_model2
gen FGPP_model3 = 1- FGP_model3

gen FGPP1 = 1- FGP1
gen FGPP2 = 1- FGP2
gen FGPP3 = 1- FGP3


sort age
#delimit ;
line FGPP_model1 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) 
|| line FGPP_model2 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(navy) 
|| line FGPP_model3 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(green) 
|| line FGPP1 age , lwidth(medthick) lcolor(cranberry)  lwidth(medthick) lpattern(dash)
|| line FGPP2  age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(navy) lpattern(dash)
|| line FGPP3  age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(green) lpattern(dash)
  ||, 
  ylabel(0.0(.03).1)
  xlabel(25(5)74)
  ytitle("Fraction")
xtitle("Age")
   title("H Transitions from Fair/Good to Poor", color(black))
  legend(label(1 "HS or Less, Model") label(2 "Some College, Model") label(3 "College, Model")  label(4 "HS or Less, Data") label(5 "Some College, Data") label(6 "College, Data"))
  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr	
graph export "${out_figures}\Health_Trans_Profile_FGP.eps", replace
}

* Figure 7
* H by age 
qui{
clear
use "Data files/Benchmark.dta"
		gen indic_P=0
		gen indic_A=0
		gen indic_G=0
		replace indic_P=1 if H==1 | H==4
		replace indic_A=1 if H==2
		replace indic_G=1 if H==3
		
replace R_recorded=R_recorded-1
replace R=R-1
		
		
		sort education age_group
		collapse (mean)   dp dp_record du du_record s s_record R R_recorded indic_P indic_A indic_G, by(education age_group) fast
        replace dp=dp*100
		replace dp_record=dp_record*100
		
		        replace du=du*100
		replace du_record=du_record*100
		
		        replace s=s*100
		replace s_record=s_record*100

		replace R_recorded=R_recorded*100
		replace R=R*100
		
		replace indic_P=indic_P*100
		replace indic_A=indic_A*100
		replace indic_G=indic_G*100
		
	
		reshape wide   dp dp_record du du_record s s_record R R_recorded indic_P indic_A indic_G, i(age_group) j(education)
		sort age_group	
		
		
merge age_group using "${mystart}\MEPS\MEPS Data Files\Intermediate Data Files\Health Moments 2 data.dta"
sort age_group
rename age_group age
*H
#delimit ;
line indic_P1 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) 
|| line  indic_A1 age , lwidth(medthick) lcolor(edkblue) lpattern(solid)
|| line  indic_G1 age , lwidth(medthick) lcolor(midgreen)
|| line H_1_data1 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) lpattern(dash)
|| line  H_2_data1 age , lwidth(medthick) lcolor(edkblue) lpattern(dash)
|| line  H_3_data1 age , lwidth(medthick) lcolor(midgreen) lpattern(dash)
  ||, 
  ylabel(0.0(10)100)
  xlabel(25(10)100)
  ytitle("Fraction")
xtitle("Age")
   title("High School or Less", color(black))
  legend(label(1 "Poor, Model") label(2 "Fair, Model") label(3 "Good, Model") label(4 "Poor, Data") label(5 "Fair, Data") label(6 "Good, Data") )
  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr
graph export "${out_figures}\Health_Profile1.eps", replace	

#delimit ;
line indic_P2 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) 
|| line  indic_A2 age , lwidth(medthick) lcolor(edkblue) lpattern(solid)
|| line  indic_G2 age , lwidth(medthick) lcolor(midgreen)
|| line H_1_data2 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) lpattern(dash)
|| line  H_2_data2 age , lwidth(medthick) lcolor(edkblue) lpattern(dash)
|| line  H_3_data2 age , lwidth(medthick) lcolor(midgreen) lpattern(dash)
  ||, 
  ylabel(0.0(10)100)
  xlabel(25(10)100)
  ytitle("Fraction")
xtitle("Age")
   title("Some College", color(black))
  legend(label(1 "Poor, Model") label(2 "Fair, Model") label(3 "Good, Model") label(4 "Poor, Data") label(5 "Fair, Data") label(6 "Good, Data") )
  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr
graph export "${out_figures}\Health_Profile2.eps", replace	


#delimit ;
line indic_P3 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) 
|| line  indic_A3 age , lwidth(medthick) lcolor(edkblue) lpattern(solid)
|| line  indic_G3 age , lwidth(medthick) lcolor(midgreen)
|| line H_1_data3 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) lpattern(dash)
|| line  H_2_data3 age , lwidth(medthick) lcolor(edkblue) lpattern(dash)
|| line  H_3_data3 age , lwidth(medthick) lcolor(midgreen) lpattern(dash)
  ||, 
  ylabel(0.0(10)100)
  xlabel(25(10)100)
  ytitle("Fraction")
xtitle("Age")
   title("College", color(black))
  legend(label(1 "Poor, Model") label(2 "Fair, Model") label(3 "Good, Model") label(4 "Poor, Data") label(5 "Fair, Data") label(6 "Good, Data") )
  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr
graph export "${out_figures}\Health_Profile3.eps", replace	

}

/*
* Figure 9 
qui{
clear
cap log close
infile using "dictionary_H.dct"
replace age=age+24

* I_treat=2 means treated. =1 not treat.
* 1=good , 2= Bad latent type        
#delimit ;
line Transition age if age<65 & H==3 & H_p==3 & education==1 & latent_health==1 & type2==2 & MC_group==2 & I_treat==1, yaxis(1) xaxis(1) lwidth(medthick)   lcolor(cranberry) lpattern(dash)
|| line Transition age if age<65 & H==3 & H_p==3 & education==1 & latent_health==2 & type2==2 & MC_group==2 & I_treat==1,  lwidth(medthick) lcolor(navy)  lpattern(dash)
|| line Transition age if age<65 & H==3 & H_p==3 & education==1 & latent_health==1 & type2==2 & MC_group==2 & I_treat==2, lwidth(medthick) lcolor(cranberry) 
|| line Transition age if age<65 & H==3 & H_p==3 & education==1 & latent_health==2 & type2==2 & MC_group==2 & I_treat==2, lwidth(medthick) lcolor(navy)
  ||, 
  ylabel(0.0(.10).90)
  xlabel(25(5)65)
  ytitle("Probability")
xtitle("Age")
   title("Probability of Keeping Good H", color(black))
   subtitle("HS or Less, Has dp shock")
  legend(label(1 "Not Treat, Good LH") label(2 "Not Treat, Bad LH") label(3 "Treat, Good LH")  label(4 "Treat, Bad LH") )
  legend(col(2) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr	
graph export "${out_figures}\H_plot1.eps", replace



#delimit ;
line Transition age if age<65 & H==3 & H_p==2 & education==1 & latent_health==1 & type2==2 & MC_group==2 & I_treat==1, yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) lpattern(dash)
|| line Transition age if age<65 & H==3 & H_p==2 & education==1 & latent_health==2 & type2==2 & MC_group==2 & I_treat==1,  lwidth(medthick) lcolor(navy)  lpattern(dash)
|| line Transition age if age<65 & H==3 & H_p==2 & education==1 & latent_health==1 & type2==2 & MC_group==2 & I_treat==2, lwidth(medthick) lcolor(cranberry) 
|| line Transition age if age<65 & H==3 & H_p==2 & education==1 & latent_health==2 & type2==2 & MC_group==2 & I_treat==2, lwidth(medthick) lcolor(navy)
  ||, 
  ylabel(0.0(.10).90)
  xlabel(25(5)65)
  ytitle("Probability")
xtitle("Age")
   title("Probability of Transitioning from Good to Fair H", color(black))
      subtitle("HS or Less, Has dp shock")
  legend(label(1 "Not Treat, Good LH") label(2 "Not Treat, Bad LH") label(3 "Treat, Good LH")  label(4 "Treat, Bad LH") )
  legend(col(2) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr	
graph export "${out_figures}\H_plot2.eps", replace



#delimit ;
line Transition age if age<65 & H==2 & H_p==1 & education==1 & latent_health==1 & type2==2 & MC_group==2 & I_treat==1, yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) lpattern(dash)
|| line Transition age if age<65 & H==2 & H_p==1 & education==1 & latent_health==2 & type2==2 & MC_group==2 & I_treat==1,  lwidth(medthick)  lcolor(navy)  lpattern(dash)
|| line Transition age if age<65 & H==2 & H_p==1 & education==1 & latent_health==1 & type2==2 & MC_group==2 & I_treat==2, lwidth(medthick) lcolor(cranberry) 
|| line Transition age if age<65 & H==2 & H_p==1 & education==1 & latent_health==2 & type2==2 & MC_group==2 & I_treat==2, lwidth(medthick)  lcolor(navy)
  ||, 
  ylabel(0.0(.10).90)
  xlabel(25(5)65)
  ytitle("Probability")
xtitle("Age")
   title("Probability of Transitioning from Fair to Poor H", color(black))
      subtitle("HS or Less, Has dp shock")
  legend(label(1 "Not Treat, Good LH") label(2 "Not Treat, Bad LH") label(3 "Treat, Good LH")  label(4 "Treat, Bad LH") )
  legend(col(2) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr	
graph export "${out_figures}\H_plot3.eps", replace
}
*/

* Figure 11
**MEDICAL CHARGES BY AGE AND INSURANCE AND EDUCATION 
qui{
clear
use "Data files/Benchmark.dta"
keep if age<65

collapse (mean)	MC, by(age_group INSURANCE education)
replace MC=MC/1000
reshape wide MC, i(age_group education) j(INSURANCE)
sort age_group education
merge age_group education using "${mystart}\MEPS\MEPS Data Files\Intermediate Data Files\Charges Moments 4 data.dta"

 #delimit ;
line MC1 age_group if education==1, yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) 
|| line  MC3 age_group if education==1, lwidth(medthick) lcolor(edkblue) 
|| line  Charges_plus_RX1 age_group if education==1, lwidth(medthick) lcolor(cranberry)   lpattern(dash)
|| line  Charges_plus_RX3 age_group if education==1, lwidth(medthick) lcolor(edkblue)  lpattern(dash)
  ||, 
  ylabel(0.0(1)12)
  xlabel(25(5)65)
  ytitle("Charges (thousands)")
xtitle("Age")
   title("Average Medical Charges, HS or Less", color(black))
  legend(label(1 "ESHI, Model") label(2 "Uninsured, Model") label(3 "ESHI, Data") label(4 "Uninsured, Data")  )
  legend(col(2) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr	
graph export "${out_figures}\Charges_profile2.eps", replace


 #delimit ;
line MC1 age_group if education==2, yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) 
|| line  MC3 age_group if education==2, lwidth(medthick) lcolor(edkblue) 
|| line  Charges_plus_RX1 age_group if education==2, lwidth(medthick) lcolor(cranberry)   lpattern(dash)
|| line  Charges_plus_RX3 age_group if education==2, lwidth(medthick) lcolor(edkblue)  lpattern(dash)
  ||, 
  ylabel(0.0(1)12)
  xlabel(25(5)65)
  ytitle("Charges (thousands)")
xtitle("Age")
   title("Average Medical Charges, Some College", color(black))
  legend(label(1 "ESHI, Model") label(2 "Uninsured, Model") label(3 "ESHI, Data") label(4 "Uninsured, Data")  )
  legend(col(2) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr	
graph export "${out_figures}\Charges_profile3.eps", replace


 #delimit ;
line MC1 age_group if education==3, yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) 
|| line  MC3 age_group if education==3, lwidth(medthick) lcolor(edkblue) 
|| line  Charges_plus_RX1 age_group if education==3, lwidth(medthick) lcolor(cranberry)   lpattern(dash)
|| line  Charges_plus_RX3 age_group if education==3, lwidth(medthick) lcolor(edkblue)  lpattern(dash)
  ||, 
  ylabel(0.0(1)12)
  xlabel(25(5)65)
  ytitle("Charges (thousands)")
xtitle("Age")
   title("Average Medical Charges, College", color(black))
  legend(label(1 "ESHI, Model") label(2 "Uninsured, Model") label(3 "ESHI, Data") label(4 "Uninsured, Data")  )
  legend(col(2) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr	
graph export "${out_figures}\Charges_profile4.eps", replace
}


* Figure 13 - Fraction not treated, and % paying bills
******* treatment and paying bills 
qui{
clear
use "Data files/Benchmark.dta"
keep if type2!=1 & age<65
		
		sort education age
		collapse (mean)   I_treat , by(education age) fast
		replace I_treat=1-I_treat
		reshape wide I_treat , i(age) j(education)
 #delimit ;
line I_treat1 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) 
|| line  I_treat2 age , lwidth(medthick) lcolor(cranberry)   lpattern(dash)
|| line  I_treat3 age , lwidth(medthick) lcolor(edkblue)  lpattern(dash)
  ||, 
  ylabel(0.0(.10)1)
  xlabel(25(5)65)
  ytitle("Fraction")
xtitle("Age")
   title("Fraction Not Treated Conditional on Health Shock", color(black))
  legend(label(1 "HS or Less") label(2 "Some College") label(3 "College")  )
  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr       
graph export "${out_figures}\Treat_Profile1.eps", replace

clear
use "Data files/Benchmark.dta"
keep if type2!=1 & age<65 & I_treat==1
		
		sort education age
		collapse (mean)   I_pay , by(education age) fast
		reshape wide I_pay , i(age) j(education)
 #delimit ;
line I_pay1 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) 
|| line  I_pay2 age , lwidth(medthick) lcolor(cranberry)   lpattern(dash)
|| line  I_pay3 age , lwidth(medthick) lcolor(edkblue)  lpattern(dash)
  ||, 
  ylabel(0.0(.10)1)
  xlabel(25(5)65)
  ytitle("Fraction")
xtitle("Age")
   title("% Paying Bills Conditional on Health Shock & Treatment", color(black))
  legend(label(1 "HS or Less") label(2 "Some College") label(3 "College")  )
  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr       
graph export "${out_figures}\Pay_Profile1.eps", replace
}

   
 *Figures 14, 15, 16, left panels
* H BY EDUC AND EMPLOYMENT
qui{
clear
use "Data files/Benchmark.dta"
*keep if latent_health==2 // good types
keep if age<65 
gen H_1=0 
replace H_1=1 if H==1
gen H_2=0 
replace H_2=1 if H==2
gen H_3=0 
replace H_3=1 if H==3
collapse H_1 H_2 H_3, by(age education employed_yn)
sort age
reshape wide H_1 H_2 H_3, i(age employed_yn) j(education)
reshape wide H*, i(age) j(employed_yn)
sort age
merge age using "${mystart}\MEPS\MEPS Data Files\Intermediate Data Files\Health Moments 6 data.dta"
drop if age>64

foreach var of varlist H_110 - H_3_data31{
	replace `var'=1-`var'
}


sort age
#delimit ;
line H_310 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) 
|| line H_311 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(navy) 
|| line H_3_data10 age , lwidth(medthick) lcolor(cranberry)  lwidth(medthick) lpattern(dash)
|| line H_3_data11  age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(navy) lpattern(dash)
  ||, 
  ylabel(0.0(.10)1)
  xlabel(25(5)64)
  ytitle("Fraction")
xtitle("Age")
  legend(label(1 "Not Employed, Model") label(2 "Employed, Model") label(3 "Not Employed, Data")  label(4 "Employed, Data") )
  legend(col(2) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr	
graph export "${out_figures}\Health_Profile16_HS.eps", replace


sort age
#delimit ;
line H_320 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) 
|| line H_321 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(navy) 
|| line H_3_data20 age , lwidth(medthick) lcolor(cranberry)  lwidth(medthick) lpattern(dash)
|| line H_3_data21  age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(navy) lpattern(dash)
  ||, 
  ylabel(0.0(.10)1)
  xlabel(25(5)64)
  ytitle("Fraction")
xtitle("Age")
  legend(label(1 "Not Employed, Model") label(2 "Employed, Model") label(3 "Not Employed, Data")  label(4 "Employed, Data") )
  legend(col(2) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr	
graph export "${out_figures}\Health_Profile16_SC.eps", replace


sort age
#delimit ;
line H_330 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) 
|| line H_331 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(navy) 
|| line H_3_data30 age , lwidth(medthick) lcolor(cranberry)  lwidth(medthick) lpattern(dash)
|| line H_3_data31  age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(navy) lpattern(dash)
  ||, 
  ylabel(0.0(.10)1)
  xlabel(25(5)64)
  ytitle("Fraction")
xtitle("Age")
  legend(label(1 "Not Employed, Model") label(2 "Employed, Model") label(3 "Not Employed, Data")  label(4 "Employed, Data") )
  legend(col(2) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr	
graph export "${out_figures}\Health_Profile16_C.eps", replace
}



*Figures 14, 15, 16, right panels
* H BY EDUC AND INSURANCE
* looking at those in fair/poor
qui{
clear
use "Data files/Benchmark.dta"
keep if age<65 
gen H_1=0 
replace H_1=1 if H==1
gen H_2=0 
replace H_2=1 if H==2
gen H_3=0 
replace H_3=1 if H==3
collapse H_1 H_2 H_3, by(age education ESHI)
sort age
reshape wide H_1 H_2 H_3, i(age ESHI) j(education)
reshape wide H*, i(age) j(ESHI)
sort age
merge age using "${mystart}\MEPS\MEPS Data Files\Intermediate Data Files\Health Moments 3 data.dta"
drop if age>64

foreach var of varlist H_110 - H_3_data31{
	replace `var' = 1-`var'
}

sort age
#delimit ;
line H_310 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) 
|| line H_311 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(navy) 
|| line H_3_data10 age , lwidth(medthick) lcolor(cranberry)  lwidth(medthick) lpattern(dash)
|| line H_3_data11  age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(navy) lpattern(dash)
  ||, 
  ylabel(0.0(.10)1)
  xlabel(25(5)64)
  ytitle("Fraction")
xtitle("Age")
  legend(label(1 "No ESHI, Model") label(2 "ESHI, Model") label(3 "No ESHI, Data")  label(4 "ESHI, Data") )
  legend(col(2) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr	
graph export "${out_figures}\Health_Profile4_HS.eps", replace


#delimit ;
line H_320 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) 
|| line H_321 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(navy) 
|| line H_3_data20 age , lwidth(medthick) lcolor(cranberry)  lwidth(medthick)  lpattern(dash)
|| line H_3_data21  age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(navy) lpattern(dash)
  ||, 
  ylabel(0.0(.10)1)
  xlabel(25(5)64)
  ytitle("Fraction")
xtitle("Age")
  legend(label(1 "No ESHI, Model") label(2 "ESHI, Model") label(3 "No ESHI, Data")  label(4 "ESHI, Data") )
  legend(col(2) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr	
graph export "${out_figures}\Health_Profile4_SC.eps", replace


#delimit ;
line H_330 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) 
|| line H_331 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(navy) 
|| line H_3_data30 age , lwidth(medthick) lcolor(cranberry)  lwidth(medthick) lpattern(dash)
|| line H_3_data31  age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(navy) lpattern(dash)
  ||, 
  ylabel(0.0(.10)1)
  xlabel(25(5)64)
  ytitle("Fraction")
xtitle("Age")
  legend(label(1 "No ESHI, Model") label(2 "ESHI, Model") label(3 "No ESHI, Data")  label(4 "ESHI, Data") )
  legend(col(2) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr	
graph export "${out_figures}\Health_Profile4_C.eps", replace
}


* Figure 17
* H by latent health type
qui{
* Here we look at H for all
clear
use "Data files/Benchmark.dta"
*keep if latent_health==2 // good types
keep if age<65 
gen H_1=0 
replace H_1=1 if H==1
gen H_2=0 
replace H_2=1 if H==2
gen H_3=0 
replace H_3=1 if H==3


collapse H_1 H_2 H_3, by(age education latent_health)
sort age
reshape wide H_1 H_2 H_3, i(age latent_health) j(education)
reshape wide H*, i(age) j(latent_health)
sort age

sort age
#delimit ;
line H_312 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) lpattern(dash)
|| line H_311 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(navy) 
  ||, 
  ylabel(0.0(.10)1)
  xlabel(25(5)64)
  ytitle("Fraction")
xtitle("Age")
   title("Fraction in Good H, HS or Less, by Latent Health Type", color(black))
  legend(label(1 "Bad Latent Health") label(2 "Good Latent Health")   )
  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr	
graph export "${out_figures}\Health_Profile13_1.eps", replace



sort age
#delimit ;
line H_322 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) lpattern(dash)
|| line H_321 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(navy)
  ||, 
  ylabel(0.0(.10)1)
  xlabel(25(5)64)
  ytitle("Fraction")
xtitle("Age")
   title("Fraction in Good H, Some College, by Latent Health Type", color(black))
  legend(label(1 "Bad Latent Health") label(2 "Good Latent Health")   )
  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr	
graph export "${out_figures}\Health_Profile14_2.eps", replace



sort age
#delimit ;
line H_332 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) lpattern(dash)
|| line H_331 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(navy)
  ||, 
  ylabel(0.0(.10)1)
  xlabel(25(5)64)
  ytitle("Fraction")
xtitle("Age")
   title("Fraction in Good H, College,  by Latent Health Type", color(black))
  legend(label(1 "Bad Latent Health") label(2 "Good Latent Health")   )
  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr	
graph export "${out_figures}\Health_Profile15_3.eps", replace
}

* Figure 18 
* ESHI by latent skill and latent health
qui{
clear
use "Data files/Benchmark.dta"
keep if age<65
collapse (mean) ESHI, by(education age productivity_type)
reshape wide ESHI, i(age education) j(productivity_type)

 #delimit ;
line ESHI1 age if education==1, yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) 
|| line  ESHI2 age if education==1, lwidth(medthick) lcolor(cranberry)   lpattern(dash)
|| line  ESHI3 age if education==1, lwidth(medthick) lcolor(edkblue)  lpattern(dash)
  ||, 
  ylabel(0.0(.10)1)
  xlabel(25(5)65)
  ytitle("Fraction")
xtitle("Age")
   title("Fraction with ESHI by Latent Skill, HS or Less", color(black))
  legend(label(1 "Low") label(2 "Medium") label(3 "High")  )
  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr       
graph export "${out_tables}\ESHI_skill1.eps", replace


 #delimit ;
line ESHI1 age if education==2, yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) 
|| line  ESHI2 age if education==2, lwidth(medthick) lcolor(cranberry)   lpattern(dash)
|| line  ESHI3 age if education==2, lwidth(medthick) lcolor(edkblue)  lpattern(dash)
  ||, 
  ylabel(0.0(.10)1)
  xlabel(25(5)65)
  ytitle("Fraction")
xtitle("Age")
   title("Fraction with ESHI by Latent Skill, Some College", color(black))
  legend(label(1 "Low") label(2 "Medium") label(3 "High")  )
  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr       
graph export "${out_tables}\ESHI_skill2.eps", replace


 #delimit ;
line ESHI1 age if education==3, yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) 
|| line  ESHI2 age if education==3, lwidth(medthick) lcolor(cranberry)   lpattern(dash)
|| line  ESHI3 age if education==3, lwidth(medthick) lcolor(edkblue)  lpattern(dash)
  ||, 
  ylabel(0.0(.10)1)
  xlabel(25(5)65)
  ytitle("Fraction")
xtitle("Age")
   title("Fraction with ESHI by Latent Skill, College", color(black))
  legend(label(1 "Low") label(2 "Medium") label(3 "High")  )
  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr       
graph export "${out_tables}\ESHI_skill3.eps", replace


*****************************************************************
*** FRACTIONs WITH ESHI BY LATENT health TYPE
*****************************************************************
clear
use "Data files/Benchmark.dta"
keep if age<65
collapse (mean) ESHI, by(education age latent_health)
reshape wide ESHI, i(age education) j(latent_health)

 #delimit ;
line ESHI2 age if education==1, yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) lpattern(dash)
|| line  ESHI1 age if education==1, lwidth(medthick) lcolor(edkblue)  
  ||, 
  ylabel(0.0(.10)1)
  xlabel(25(5)65)
  ytitle("Fraction")
xtitle("Age")
   title("Fraction with ESHI by Latent Health, HS or Less", color(black))
  legend(label(1 "Bad Latent Health") label(2 "Good Latent Health")   )
  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr       
graph export "${out_tables}\ESHI_LH1.eps", replace



 #delimit ;
line ESHI2 age if education==2, yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) lpattern(dash)
|| line  ESHI1 age if education==2, lwidth(medthick) lcolor(edkblue) 
  ||, 
  ylabel(0.0(.10)1)
  xlabel(25(5)65)
  ytitle("Fraction")
xtitle("Age")
   title("Fraction with ESHI by Latent Health, Some College", color(black))
  legend(label(1 "Bad Latent Health") label(2 "Good Latent Health")   )
  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr       
graph export "${out_tables}\ESHI_LH2.eps", replace


 #delimit ;
line ESHI2 age if education==3, yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) lpattern(dash)
|| line  ESHI1 age if education==3, lwidth(medthick) lcolor(edkblue) 
  ||, 
  ylabel(0.0(.10)1)
  xlabel(25(5)65)
  ytitle("Fraction")
xtitle("Age")
   title("Fraction with ESHI by Latent Health, College", color(black))
  legend(label(1 "Bad Latent Health") label(2 "Good Latent Health")   )
  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr       
graph export "${out_tables}\ESHI_LH3.eps", replace
}


* Figure 20	
*EMP FT by H
qui{
clear
use "Data files/Benchmark.dta"
		keep if age<65
replace emp_ft_pt=0 if emp_ft_pt==1 
replace emp_ft_pt=1 if emp_ft_pt==2 // so now the FT == 1 and all else ==0
	
collapse (mean) emp_ft_pt, by(education age H)
	
rename education EDUCTYPE
rename age AGE
rename emp_ft_pt number_model
replace number_model=0 if number_model==.
replace number_model= number_model*100

reshape wide number_model, i(AGE EDUCTYPE) j(H)
rename AGE age
sort age
reshape wide number*, i(age) j(EDUCTYPE)
sort age


merge age using "${mystart}\MEPS\MEPS Data Files\Intermediate Data Files\Employment Distribution Health No Shocks.dta"
sort age
*we only want the FT
		#delimit ;
		line Employment11 age, yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) lpattern(dot)
		|| line Employment21 age , lwidth(medthick) lcolor(midgreen) lpattern(shortdash)
		|| line Employment31 age , lwidth(medthick) lcolor(edkblue) lpattern(shortdash)
		|| line  number_model11 age , lwidth(medthick)  lcolor(cranberry)
		|| line  number_model21 age , lwidth(medthick) lcolor(midgreen)
		|| line  number_model31 age , lwidth(medthick) lcolor(edkblue)
		  ||, 
		  ylabel(0.0(20)100)
		  xlabel(25(10)64)
		  ytitle("Percentage")
		   xtitle("Age")
		  title("High School or Less", color(black))
		  legend(label(1 "Poor, Data") label(2 "Fair, Data") label(3 "Good, Data")  label(4 "Poor, Model") label(5 "Fair, Model") label(6 "Good, Model")       )
	  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
		  ;
		#delimit cr
graph export "${out_figures}\Emp_FT_Health1.eps", replace


		#delimit ;
		line Employment12 age, yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) lpattern(dot)
		|| line  Employment22 age , lwidth(medthick) lcolor(midgreen) lpattern(shortdash)
		|| line Employment32 age , lwidth(medthick) lcolor(edkblue) lpattern(shortdash)
		|| line  number_model12 age , lwidth(medthick)  lcolor(cranberry)
		|| line  number_model22 age , lwidth(medthick) lcolor(midgreen)
		|| line  number_model32 age , lwidth(medthick) lcolor(edkblue)
		  ||, 
		  ylabel(0.0(20)100)
		  xlabel(25(10)64)
		  ytitle("Percentage")
		   xtitle("Age")
		  title("Some College", color(black))
		  legend(label(1 "Poor, Data") label(2 "Fair, Data") label(3 "Good, Data")  label(4 "Poor, Model") label(5 "Fair, Model") label(6 "Good, Model")       )
	  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
		  ;
		#delimit cr
graph export "${out_figures}\Emp_FT_Health2.eps", replace


		#delimit ;
		line Employment13 age, yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) lpattern(dot)
		|| line  Employment23 age , lwidth(medthick) lcolor(midgreen) lpattern(shortdash)
		|| line Employment33 age , lwidth(medthick) lcolor(edkblue) lpattern(shortdash)
		|| line  number_model13 age , lwidth(medthick)  lcolor(cranberry)
		|| line  number_model23 age , lwidth(medthick) lcolor(midgreen)
		|| line  number_model33 age , lwidth(medthick) lcolor(edkblue)
		  ||, 
		  ylabel(0.0(20)100)
		  xlabel(25(10)64)
		  ytitle("Percentage")
		   xtitle("Age")
		  title("College", color(black))
		  legend(label(1 "Poor, Data") label(2 "Fair, Data") label(3 "Good, Data")  label(4 "Poor, Model") label(5 "Fair, Model") label(6 "Good, Model")       )
	  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
		  ;
		#delimit cr
graph export "${out_figures}\Emp_FT_Health3.eps", replace

}


* Figure 21 - Emp by PT and PT
* Figure 22 - Emp by educ and hours and insurance
* Figure 23 - emp in CPS and MEPS
* EMPLOYMENT
qui{
clear
use "Data files/Benchmark.dta"
keep if age<65

tab employment_offer_yn employed_yn
 
replace emp=2 if emp==3 // PT
replace emp=3 if emp==4 // FT, no ins
replace emp=4 if emp==5 // FT, ins
collapse  (count) number_model=ID, by(emp education age)
reshape wide number_model, i(age education) j(emp)
egen total=rowtotal(number_model1 number_model2 number_model3 number_model4)
replace number_model1=number_model1/total*100
replace number_model2=number_model2/total*100
replace number_model3=number_model3/total*100
replace number_model4=number_model4/total*100
sort age education
merge age education using  "${mystart}\CPS\CPS Data Files\Intermediate Data Files\CPS LS Ins.dta"
drop _merge
foreach num of numlist 1 2 3 4 {
    rename number`num' number`num'_CPS
}
drop number5 // this is the same as 4
cap drop total
sort age education

merge age education using "${mystart}\MEPS\MEPS Data Files\Intermediate Data Files\MEPS LS Ins.dta"


sort age
#delimit ;
line number2_CPS age if education==1, yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) lpattern(dash)
|| line  number3_CPS age if education==1, lwidth(medthick) lcolor(edkblue) lpattern(dash)
|| line  number4_CPS age if education==1, lwidth(medthick) lcolor(green) lpattern(dash)
|| line  number2 age if education==1, lwidth(medthick) lcolor(cranberry)  
|| line  number3 age if education==1, lwidth(medthick) lcolor(edkblue) 
|| line  number4 age if education==1, lwidth(medthick) lcolor(green) 
  ||, 
  ylabel(0.0(10)90)
  xlabel(25(10)65)
  ytitle("Percentage")
xtitle("Age")
  title("Distribution of Employment", color(black))
subtitle("High School or Less")
  legend(label(1 "PT, CPS") label(2 "FT, no ESHI, CPS") label(3 "FT, ESHI, CPS")   label(4 "PT, MEPS") label(5 "FT, no ESHI, MEPS") label(6 "FT, ESHI, MEPS") )
  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr
graph export "${out_figures}\LF1_CPS_MEPS.eps", replace


sort age
#delimit ;
line number2_CPS age if education==2, yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry)  lpattern(dash)
|| line  number3_CPS age if education==2, lwidth(medthick) lcolor(edkblue)  lpattern(dash)
|| line  number4_CPS age if education==2, lwidth(medthick) lcolor(green)  lpattern(dash)
|| line  number2 age if education==2, lwidth(medthick) lcolor(cranberry)  
|| line  number3 age if education==2, lwidth(medthick) lcolor(edkblue) 
|| line  number4 age if education==2, lwidth(medthick) lcolor(green) 
  ||, 
  ylabel(0.0(10)90)
  xlabel(25(10)65)
  ytitle("Percentage")
xtitle("Age")
  title("Distribution of Employment", color(black))
subtitle("Some College")
  legend(label(1 "PT, CPS") label(2 "FT, no ESHI, CPS") label(3 "FT, ESHI, CPS")   label(4 "PT, MEPS") label(5 "FT, no ESHI, MEPS") label(6 "FT, ESHI, MEPS") )
  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr
graph export "${out_figures}\LF2_CPS_MEPS.eps", replace



sort age
#delimit ;
line number2_CPS age if education==3, yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) lpattern(dash)
|| line  number3_CPS age if education==3, lwidth(medthick) lcolor(edkblue) lpattern(dash)
|| line  number4_CPS age if education==3, lwidth(medthick) lcolor(green) lpattern(dash)
|| line  number2 age if education==3, lwidth(medthick) lcolor(cranberry) 
|| line  number3 age if education==3, lwidth(medthick) lcolor(edkblue) 
|| line  number4 age if education==3, lwidth(medthick) lcolor(green) 
  ||, 
  ylabel(0.0(10)90)
  xlabel(25(10)65)
  ytitle("Percentage")
xtitle("Age")
  title("Distribution of Employment", color(black))
subtitle("College")
  legend(label(1 "PT, CPS") label(2 "FT, no ESHI, CPS") label(3 "FT, ESHI, CPS")   label(4 "PT, MEPS") label(5 "FT, no ESHI, MEPS") label(6 "FT, ESHI, MEPS") )
  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr
graph export "${out_figures}\LF3_CPS_MEPS.eps", replace



sort age
#delimit ;
line number_model2 age if education==1, yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry)
|| line  number_model3 age if education==1, lwidth(medthick) lcolor(edkblue) 
|| line  number_model4 age if education==1, lwidth(medthick) lcolor(green) 
|| line  number2_CPS age if education==1, lwidth(medthick) lcolor(cranberry)  lpattern(dash)
|| line  number3_CPS age if education==1, lwidth(medthick) lcolor(edkblue) lpattern(dash)
|| line  number4_CPS age if education==1, lwidth(medthick) lcolor(green) lpattern(dash)
  ||, 
  ylabel(0.0(10)90)
  xlabel(25(10)65)
  ytitle("Percentage")
xtitle("Age")
  title("Distribution of Employment", color(black))
subtitle("High School or Less")
  legend(label(1 "PT, Model") label(2 "FT, no ESHI, Model") label(3 "FT, ESHI, Model")   label(4 "PT, Data") label(5 "FT, no ESHI, Data") label(6 "FT, ESHI, Data") )
  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr
graph export "${out_figures}\LF1_ESHI.eps", replace

sort age
#delimit ;
line number_model2 age if education==2, yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry)
|| line  number_model3 age if education==2, lwidth(medthick) lcolor(edkblue) 
|| line  number_model4 age if education==2, lwidth(medthick) lcolor(green) 
|| line  number2_CPS age if education==2, lwidth(medthick) lcolor(cranberry)  lpattern(dash)
|| line  number3_CPS age if education==2, lwidth(medthick) lcolor(edkblue) lpattern(dash)
|| line  number4_CPS age if education==2, lwidth(medthick) lcolor(green) lpattern(dash)
  ||, 
  ylabel(0.0(10)90)
  xlabel(25(10)65)
  ytitle("Percentage")
xtitle("Age")
  title("Distribution of Employment", color(black))
subtitle("Some College")
  legend(label(1 "PT, Model") label(2 "FT, no ESHI, Model") label(3 "FT, ESHI, Model")   label(4 "PT, Data") label(5 "FT, no ESHI, Data") label(6 "FT, ESHI, Data") )
  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr
graph export "${out_figures}\LF2_ESHI.eps", replace


sort age
#delimit ;
line number_model2 age if education==3, yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry)
|| line  number_model3 age if education==3, lwidth(medthick) lcolor(edkblue) 
|| line  number_model4 age if education==3, lwidth(medthick) lcolor(green) 
|| line  number2_CPS age if education==3, lwidth(medthick) lcolor(cranberry)  lpattern(dash)
|| line  number3_CPS age if education==3, lwidth(medthick) lcolor(edkblue) lpattern(dash)
|| line  number4_CPS age if education==3, lwidth(medthick) lcolor(green) lpattern(dash)
  ||, 
  ylabel(0.0(10)90)
  xlabel(25(10)65)
  ytitle("Percentage")
xtitle("Age")
  title("Distribution of Employment", color(black))
subtitle("College")
  legend(label(1 "PT, Model") label(2 "FT, no ESHI, Model") label(3 "FT, ESHI, Model")   label(4 "PT, Data") label(5 "FT, no ESHI, Data") label(6 "FT, ESHI, Data") )
  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr
graph export "${out_figures}\LF3_ESHI.eps", replace



* NOW LET'S COMBINE THOSE WITH AND WITHOUT INSURANCE WORKING FT
replace number_model3 = number_model3+ number_model4
drop number_model4
replace number3 = number3+ number4
drop number4
replace number3_CPS = number3_CPS+ number4_CPS
drop number4_CPS

sort age
#delimit ;
line number_model2 age if education==1, yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry)
|| line  number_model3 age if education==1, lwidth(medthick) lcolor(edkblue) 
|| line  number2_CPS age if education==1, lwidth(medthick) lcolor(cranberry)  lpattern(dash)
|| line  number3_CPS age if education==1, lwidth(medthick) lcolor(edkblue) lpattern(dash)
  ||, 
  ylabel(0.0(10)90)
  xlabel(25(10)65)
  ytitle("Percentage")
xtitle("Age")
  title("High School or Less", color(black))
  legend(label(1 "PT, Model") label(2 "FT, Model")   label(3 "PT, Data") label(4 "FT, Data")  )
  legend(col(2) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr
graph export "${out_figures}\LF1.eps", replace

sort age
#delimit ;
line number_model2 age if education==2, yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry)
|| line  number_model3 age if education==2, lwidth(medthick) lcolor(edkblue) 
|| line  number2_CPS age if education==2, lwidth(medthick) lcolor(cranberry)  lpattern(dash)
|| line  number3_CPS age if education==2, lwidth(medthick) lcolor(edkblue) lpattern(dash)
  ||, 
  ylabel(0.0(10)90)
  xlabel(25(10)65)
  ytitle("Percentage")
xtitle("Age")
  title("Some College", color(black))
  legend(label(1 "PT, Model") label(2 "FT, Model")   label(3 "PT, Data") label(4 "FT, Data")  )
  legend(col(2) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr
graph export "${out_figures}\LF2.eps", replace


sort age
#delimit ;
line number_model2 age if education==3, yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry)
|| line  number_model3 age if education==3, lwidth(medthick) lcolor(edkblue) 
|| line  number2_CPS age if education==3, lwidth(medthick) lcolor(cranberry)  lpattern(dash)
|| line  number3_CPS age if education==3, lwidth(medthick) lcolor(edkblue) lpattern(dash)
  ||, 
  ylabel(0.0(10)90)
  xlabel(25(10)65)
  ytitle("Percentage")
xtitle("Age")
  title("College", color(black))
  legend(label(1 "PT, Model") label(2 "FT, Model")   label(3 "PT, Data") label(4 "FT, Data")  )
  legend(col(2) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr
graph export "${out_figures}\LF3.eps", replace

}


* Figure 24
*VARIANCE OF EARNINGS
qui{
clear
use "Data files/Benchmark.dta"

drop if age>64
*keep if H==3
replace wage_noise=. if wage_noise==0
gen wages=wage_noise
replace wages=94.34 if wages>94.34 & wages!=. // to make the same as MEPS where wages top coded
		*gen emp_type=. 
		*replace emp_type=1 if (offer_accepted==2 | offer_accepted==3) // PT
		*replace emp_type=2 if (offer_accepted==4 | offer_accepted==5) // FT
		keep if offer_accepted==4 | offer_accepted==5
		gen lwages=log(wages)
keep if H==3 // as in MEPS		
		collapse (p10) p10=lwages (p90) p90=lwages  (sd) var_wages =lwages, by(education age) 
		replace var_wages=var_wages^2
		rename var_wages var_wages_model
		gen ratio_90_10 = p90-p10
		drop p10 p90
		reshape wide var_wages_model ratio_90_10, i(age) j(education)
		
		sort age
		merge age using "${mystart}\MEPS\MEPS Data Files\Intermediate Data Files\Wages Good H.dta"
		drop _merge
		
		sort age
		drop if age<25
		
		#delimit ;
		line var_wages_model1 age, yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) lpattern(dot)
		|| line  var_wages1 age , lwidth(medthick) lcolor(edkblue)
		||, 
				  ylabel(0.0(.1).5)
				  xlabel(25(5)64)
				  ytitle("Variance of Log Wages, FT Workers")
				  xtitle("Age")
				  title("Variance of Log Wages, FT Workers", color(black))
				subtitle("High School or Less")
  legend(label(1 "Model") label(2 "Data") )
  legend(col(2) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
				  ;
				#delimit cr
graph export "${out_figures}\Var_Wages1.eps", replace

		#delimit ;
		line var_wages_model2 age, yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) lpattern(dot)
		|| line  var_wages2 age , lwidth(medthick) lcolor(edkblue)
		||, 
				  ylabel(0.0(.1).5)
				  xlabel(25(5)64)
				  ytitle("Variance of Log Wages, FT Workers")
				  xtitle("Age")
				  title("Variance of Log Wages, FT Workers", color(black))
				subtitle("Some College")
  legend(label(1 "Model") label(2 "Data") )
  legend(col(2) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
				  ;
				#delimit cr
graph export "${out_figures}\Var_Wages2.eps", replace


		#delimit ;
		line var_wages_model3 age, yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) lpattern(dot)
		|| line  var_wages3 age , lwidth(medthick) lcolor(edkblue)
		||, 
				  ylabel(0.0(.1).5)
				  xlabel(25(5)64)
				  ytitle("Variance of Log Wages, FT Workers")
				  xtitle("Age")
				  title("Variance of Log Wages, FT Workers", color(black))
				subtitle("College")
  legend(label(1 "Model") label(2 "Data") )
  legend(col(2) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
				  ;
				#delimit cr
graph export "${out_figures}\Var_Wages3.eps", replace
}


* Figure 25 right panel 
* Gov Tr by educ and age
qui{
clear
use "Data files/Benchmark.dta"
keep if  age<64
		
		sort education age
		collapse (mean)   TR , by(education age) fast
		reshape wide TR , i(age) j(education)
 #delimit ;
line TR1 age , yaxis(1) xaxis(1) lwidth(medthick) lcolor(cranberry) 
|| line  TR2 age , lwidth(medthick) lcolor(cranberry)   lpattern(dash)
|| line  TR3 age , lwidth(medthick) lcolor(edkblue)  lpattern(dash)
  ||, 
  ylabel(0.0(.02).2)
  xlabel(25(5)65)
  ytitle("Fraction")
xtitle("Age")
   title("Fraction getting Gov Transfers, Model", color(black))
  legend(label(1 "HS or Less") label(2 "Some College") label(3 "College")  )
  legend(col(3) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
  ;
#delimit cr       
graph export "${out_figures}\TR_Profile.eps", replace
}


* Figure 28
* earnings inequality by age
qui{
clear
use "Data files/Benchmark.dta"
keep if age<65

cap drop earnings_noise  // we did this wrong initially when reading data, now we do it right
* Eearnings with noise
gen earnings_noise =0 if offer_accepted==1
replace earnings_noise =wage_noise*hours*52 if offer_accepted!=1
rename earnings_noise incwage_CPI_adjusted
rename wage_noise incwage_hourly
rename offer_accepted labor_force



gen GINI=.
tostring GINI , replace
foreach a of numlist 25/64  {
inequal7 incwage_CPI_adjusted  if  age==`a'
replace GINI=r(gini) if  age==`a'
}
destring   GINI, replace
rename GINI GINI_model

collapse (mean)  GINI_model , by(age) 


sort age
merge age using "${mystart}\CPS\CPS Data Files\Intermediate Data Files\Gini.dta"

									#delimit ;
				 line  GINI age , lwidth(medthick) lcolor(edkblue) lpattern(dash) 
			   || line  GINI_model age  , lwidth(medthick) lcolor(cranberry)  

				  ||, 
				  ylabel(0.0(.1).8)
				  xlabel(25(5)65)
				  ytitle("Earnings Gini")
				  xtitle("Age")
				  title("Earnings Gini Coefficients, by Age", color(black) size(medium))
				  legend(label(1 "CPS Data")  label(2 "Model")   )
  legend(col(2) pos(6) region(lcolor(gs16)))
plotregion(margin(r+7 l+5) style(none))
  graphregion(icolor(white) fcolor(gs14) margin(none ))
  graphregion(color(white)) bgcolor(white)
				  ;
				#delimit cr
graph export "${out_figures}\Gini.eps", replace
}







